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1 Introduction 



00 ' We have computed the next-to-next-to-leading-order (NNLO) contributions to the 

£\J i evolution of unpolarized parton distributions in perturbative QCD [1,2]. In this 

talk, we briefly recall why this huge computation was necessary and outline how 
it was performed. We then illustrate the structure of the results and discuss their 
end-point limits which include the three-loop cusp anomalous dimensions of the 
Wilson lines. Finally the numerical impact of the new contributions is illustrated. 

(N ' 
CO 

t> 

O . 

, For the next decade, the highest-energy experiments in particle physics will be done 

f^S ' at the (anti-)proton-proton colliders Tevatron and LHC. At such machines, if 

we disregard power corrections and observables involving final-state fragmentation 
O j' functions, the cross sections for hard processes h can be schematically written as 

a: < p = £w>*f • ( L1 ) 



/,/' 



Here f p stands for the universal momentum distributions of the partons / in the 
proton, / = qi,qi,g with i = 1, . . . ,rif, where is the number of effectively massless 
. quark flavours, &[ l represent the hard (partonic) cross sections for the process 

^ ' under consideration. Hence quantitative studies of the standard model, and of 

expected and unexpected new particles, require a precise understanding of the par- 
tonic luminosities and of the QCD corrections to the corresponding cross sections. 

For many important processes, like Higgs-boson production, the second-order 
(NNLO) QCD corrections need to be taken into account, i.e., the third term in 



~(0) , -(1) , 2 -(2) 
a h + fl s °h + a s <* h 



(1.2) 



The consistent inclusion of in Eq. requires parton distributions evolved 

with the corresponding (process-independent) NNLO splitting functions 

P f ?, NLO = a s 4°) + a»P$ + alP$ . (1.3) 

The one- and two-loop splitting functions have been known for a long time [3*|-|llj. 
For the three-loop splitting functions P^ 2 \ on the other hand, only partial results 
had been obtained until recently C2] — E3 - However earlier this year we have, finally, 
computed the complete expressions of these functions |1I2| . 
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2 Outline of the calculation 



We have derived the NNLO splitting functions by computing the partonic structure 
functions in inclusive deep-inelastic scattering (DIS), 7*(g) + f(p) — > X with Q 2 = 
—q 2 > and p 2 = 0, up to the third order in the strong coupling a s = a s /(47r). This 
computation has been performed for all even or odd values of the Mellin variable 
N via the three-loop forward Compton amplitudes. -f*(q) + f(p) — > J*(q) + f(p)- 

This approach has two major advantages: Firstly it enables us to obtain, at 
almost the same time, also the three-loop coefficient functions in DIS [23! ■ Secondly 
it allows us to check our programs, at almost any stage, by falling back to the 
Mincer program |24I25| employed in the fixed- N calculations of refs. [r^-pH]. 



2.1 Mass-factorization in DIS 

Before we address the main computational task, we briefly sketch how the splitting 
functions are extracted from the calculation. We start by writing the physical 
structure functions F a in terms of the (perturbatively calculable) bare partonic 
structure functions F aik , the bare coupling a s and the bare parton distributions f k , 

F a (Q 2 ) = F a , k (a s ,Q 2 ,e)*f k . (2.1) 

Summation over the parton species k is understood, and * stands for either the 
convolution in Bjorkcn-a; space or a simple multiplication of the Mellin moments. 
As indicated in Eq. (|2.1|l . we use dimensional regularization with D = 4 — 2s, thus 
the singularities of F a k appear as poles e~ l . After the ultraviolet divergences have 
been removed by coupling-constant renormalization, at the renormalization scale 
fi r , only initial-state mass singularities remain. They arise when two momenta 
become collinear, e.g., p and k in Fig. 1, leading to propagator denominators 

(p-k) 2 = -2\p\\k\(l-cos$) -\p\\k\$ 2 . 

These singularities are removed by mass factorization, at the factorization scale fi / : 
F a ,k is decomposed into finite pieces, the coefficient functions C a ,i, and the universal 
transition functions which contain the (a-independent) pole parts of F a k . The 
latter are combined with the f k to form the finite renormalized parton densities fi, 

F a (Q 2 ) = F aik (a s (ti 2 r ),%,e)*f k fi(j* 2 f ,l4) 

Mr a 



The decomposition (|2.2|l is not unique. We employ the usual MS scheme where, 
besides the 1/e poles, only the artefacts S e = exp(e[ln(47r)— 7 e ]) of dimensional reg- 
ularization are removed from the coefficient functions. Differentiation of Eq. I|2.2|) 
finally leads to the evolution equations for the renormalized parton distributions, 

-J-^U - Jr^ */* = #\ *^*f 3 = P*i*f 3 ■ (2-3) 

olnjiif dm/if ohifXf J 

The splitting functions (|1.3(l can thus be obtained from the 1/e poles in Eq. Q2.2[l . 



C a , l (a s (fi 2 r ),^,^)*T lk (a s ( f i 2 r ),^,e)*f k . (2.2) 
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Figure 1. Sample diagrams contributing to the inclusive process -y*g — * X up to third order in a s 



2.2 The flavour decomposition 

It is convenient to decompose the system <|2.3[) of 2nt + 1 coupled equation as far 
as possible from charge conjugation and flavour symmetry constraints alone. The 
general structure of the (anti-)quark (anti-)quark splitting functions reads 



Pa 



Pn 



■Sr,P v + P s 
- "ifc^qq T 1 qq 



■■s ik p: 



p* 



(2.4) 



This structure leads to three independently evolving types of flavour non-singlet 
combinations. The flavour asymmetries q^j, and the total valence distribution q^ 



S-ns.ik 



S, ± ?i - ilk ± Sk) , 9ns = T,rU iSr 

respectively evolve with 



P 



pv ± pv 

qq qq 



pv _ pv _ pv 

ns qq qq 



The singlet quark distribution, q, 

d fq h 
S 



n fiPqq Pqq) ~ 



P E 



ns ' 

(2.5) 



(2.6) 



Sr=i (9r + 9r) ^ s coupled to gluon density g, 

(2-7) 



P P 
-'qq ^qg 

P P 

^gq ^gg 



din /Xy 

where the quark-quark splitting function P qq can be expressed as 



Pqq — P ns + n fiPqq + -Pqq) — ^ns + P\ 



The off-diagonal entries in Eq. I|2.7|l are given by 



-Pqg — P qig , 



P — P 

-*gq ~~ -^gq; 



(2.8) 



(2.9) 

p - 



in terms of the flavour-independent splitting functions P qig = Pq ig and P gti 

In the expansion in powers of a s , the flavour-diagonal ('valence') quantity P q q in 
Eq. 12. 4[) starts at first order. P q ^ and the flavour-independent ('sea') contributions 
P qq and P qq - and hence the 'pure-singlet' term P ps - are of order a?. A non- 



vanishing P n s s = P qq — P qq in Eq. I|2.6|) occurs for the first time at the third order 
and introduces a new colour structure, d a bcd a bc- See Fig. lb of ref. [2U for a typical 
diagram contributing to P n s s . 
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2.3 Set-up of the calculation 

Diagrams like those in Fig. 1 have been calculated directly, by working out the 
phase-space integrations, for the derivation of the complete second-order coefficient 
functions EZj-[20]. An extension of this procedure to the third order, however, 
does not seem feasible. Instead, we employ the optical theorem 




to transform the problem into forward Compton amplitudes. We then make use of 
a theorem that the coefficient of (2p ■ q) N provides the iV-th Mellin moment, 

F(N) = [ dxx N - 1 F(x) , 
Jo 

of the partonic structure functions (|2.1fl which we need to calculate. 

In order to obtain the complete set of the third-order contributions to the split- 
ting functions (|2.6|) and 12. 7|) we have to include, besides the photon shown above, 
also the IV-boson [25] - for accessing P n ~ and P n s s - and a fictitious classical scalar 
<f) coupling directly only to the gluon field via 4>G^ v G^ y |32I33| - for accessing 
P gq and P gg . Especially the latter leads to a substantial increase of the number of 
diagrams (generated with Qgraf j^l]) as shown in Table 1. Among the partons / 
we also include the standard ghost h. This allows us to take the sum over external 
gluon spins by contracting with — g M „ instead of the full physical expression. 





tree 


1-loop 


2-loop 


3-loop 


<Y) 


1 


3 


25 


359 


gl 




2 


17 


345 


h-j 






2 


56 


qW 


1 


3 


32 


589 


q<* 




1 


23 


696 




1 


8 


218 


6378 


h(f> 




1 


33 


1184 


sum 


3 


18 


350 


9607 



Table 1. The number of diagrams employed in our calculation of the three-loop splitting functions. 

Obviously a highly efficient symbolic treatment is needed to cope with the task at 
hand. Unsurprisingly, we use FORM for all manipulations. Note that the capabili- 
ties of this program had to be extended substantially for this computation |35I36| . 
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2-4 Treatment of the integrals 

Finally we illustrate the computation of the integrals required to evaluate the for- 
ward Compton amplitudes. One of the 9607 three-loop diagrams in Table 1 is 
shown here together with a useful pictorial representation of its momentum flow: 




For the latter we temporarily disregard the external parton lines and draw the 
remaining self-energy type diagram, the topology of which is denoted following the 
notation of ref. |25J ■ Our example is a ladder (LA) diagram. The (partly additional) 
denominators carrying the parton momentum p are then indicated by the fat (in the 
coloured version: red) lines. Here p runs, after turning the diagram upside-down, 
through the lines 1, 2 and 3, thus the example is assigned the subtopology LA13. 

According to our discussion in the previous subsection, we need analytic expres- 
sions for the (dimensionless) coefficients I(N) of (2p ■ q) /Q 2a - One might try to 
obtain I(N) by brute force, Taylor-expanding the denominators withp and working 
out the sums. It turns out that such a strategy, in general, does not work. Instead, 
we employ identities based on integration by parts, scaling arguments and form- 
factor decompositions (see Sect. 2 of ref. pQ) to successively simplify the integrals. 

The LA13 integrals, e.g., can be simplified by applying p^d/dq^ both inside and 
outside the integral. For the scalar integral with unit denominators this yields 




Here the LA13 integral occurs twice, once with a prefactor 2p ■ q. Hence Eq. (2.10) 
represents a difference equation (here of order n = 1) which expresses its coefficient 
I{N) in terms of that of a LA12 integral with an enhanced denominator in the 3- line, 

a (N)I{N) - ... - a n (N)I(N-n) - G(N) = . (2.11) 

First-order recursion relations like Eq. (2.10) can be reduced to a sum. Higher- 
order recursions (we have used equations up to n = 4) can be solved by inserting 
a suitable ansatz into Eq. (|2.11|l . Both procedures make use of the fact that all 
integrals required for the computation of the splitting functions can be expressed 
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in terms of harmonic sums |37| . Recall that these sums are recursively defined by 
S± m (M) = ' S ±mum2 ,..., mk (M) = £ ^5 m „..., mfc (i) . (2.12) 

i=l i=l 

To the accuracy in the dimensional offset e required for the calculation of the 
splitting functions, our example integral reads, using the FORM notations den(i+N) 
for l/(N + i) and S(R(ml, . . . ,mk) ,i+N) for S mu .... mk (N+i), 

l 1 1 



(2.13) 



+theta(N)*sign(N)*ep--2*(-8/3*den(l+N)-2-4*den(l+N)~3+8/3*den(l+N)~2*S(R( 

1) ,l+N)+4/3*den(i+N)*S(R(l) , l+N)+2/3*den(l+N) *S(R(2) , l+N)-4/3*den(2+N) 
~2-2*den(2+N)~3+4/3*den(2+N)~2*S(R(l) ,2+N)+4/3*den(2+N)*S(R(l) ,2+N)+2/ 
3*den(2+N)*S(R(2) ,2+N)+4/3*S(R(l) ,N)+2/3*S(R(l , 2) ,N)-2*S(R(2) ,N)-4/3* 
S(R(2) ,N)*N+4*S(R(2,1) ,N)+4/3*S(R(2, 1) ,N)*N-6*S(R(3) ,N)-2*S(R(3) ,N)*N) 

+theta(N) *sign(N) *ep~-l* (32*den( 1+N) ~2+164/3*den(l+N) ~3+24*den( 1+N) "4-20/ 
3*den(l+N)~3*S(R(l) ,l+N)-88/3*den(l+N)~2*S(R(l) , l+N)+8/3*den(l+N) ~2*S ( 
R(l,l) , l+N)-40/3*den(l+N)~2*S(R(2) , l+N)-16*den(l+N) *S(R(1) , 1+N) +8/3* 
den(l+N)*S(R(l,l) ,l+N)+10/3*den(l+N)*S(R(l,2) ,l+N)-58/3*den(l+N)*S(R(2 
) ,l+N)+10*den(l+N)*S(R(2,l) ,l+N)-18*den(l+N)*S(R(3) ,l+N)+16*den(2+N)~2 
+82/3*den(2+N)~3+12*den(2+N)~4-10/3*den(2+N)~3*S(R(l) , 2+N) -44/3*den(2+ 
N)"2*S(R(1) ,2+N)-6*den(2+N)-2*S(R(2) ,2+N)-16*den(2+N)*S(R(l) ,2+N)+8/3* 
den(2+N) *S (R(l , 1) , 2+N)+10/3*den(2+N) *S (R(l , 2) , 2+N) -46/3*den(2+N) *S (R(2 
) ,2+N)+6*den(2+N)*S(R(2,l) ,2+N)-12*den(2+N)*S(R(3) ,2+N)-20*S(R(l) ,N)+8/ 
3*S(R(1,1) ,N)+10/3*S(R(1,1,2) ,N) -16*S(R(1 , 2) ,N) -4*S (R(l , 2) ,N)*N+14*S( 
R(l,2,l) ,N)+4*S(R(1,2,1) ,N)*N-24*S(R(1,3) ,N)-6*S(R(1 ,3) , N) *N+56/3*S(R( 

2) ,N)+20*S(R(2) ,N)*N-134/3*S(R(2,1) ,N) -56/3*S(R(2 , 1) ,N)*N+16/3*S(R(2, 1 
,1) ,N)+8/3*S(R(2,l,l) ,N)*N-62/3*S(R(2,2) ,N)-22/3*S(R(2,2) ,N) *N+76*S(R( 

3) ,N)+100/3*S(R(3) ,N)*N-10*S(R(3,1) ,N) -10/3*S(R(3 , 1) ,N)*N+36*S(R(4) ,N) 
+12*S(R(4) ,N)*N) . 

Despite being uncharacteristically simple in both derivation and size, Eq. (2.10) 
illustrates the strict hierarchy of subtopologies in our procedure. Our LA13 example 
can only be evaluated once the LA12 integral in Eq. (2.10) is known. This integral, 
in turn, requires the so-called basic building blocks (with only one p-dependent 
denominator) LAn and LA22 together with other integrals of simpler topologies 
where one of the non-p denominators has been removed. Also those integrals need 
to be evaluated in terms of yet simpler cases, and so on. 

Constructing the reduction chains for all subtopologies, and computing all in- 
tegrals required for evaluating either diagrams or other, higher-level integrals took 
literally years of both human and computing resources. It would not have been pos- 
sible to get through without extensive tabulation of intermediate results for which 
new features were added to Form [201. At the end, a database had been accumu- 
lated of more than 100 000 integrals requiring about 3.5 GBytes of disk space. 
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3 Sample results in AT-space and :c-space 

We illustrate our final results by writing down the even-N anomalous dimensions 
and the corresponding x-spacc splitting functions in Quantum-Gluodynamics, i.e., 
for QCD with rif = quark flavours. The complete QCD results in refs. |1I2| are 
considerably longer, by a factor of about 15, but not structurally more complicated. 

3.1 Expressions in Mellin-N space 



We start in iV-space where, as discussed above, the actual calculations have been 
performed. Recall that we expand in terms of a s = a s /(47r), and that ■y^(N) = 
! — p( n \N) is the coefficient of a" +1 . Here we hide ^-dependent denominators 
by using differences of harmonic sums at suitably shifted arguments for which we 
employ the abbreviations 

N±iS A = S A (N±i) , N±=N ±1 . 

In this notation the well-known one- and two-loop results |3l4l7ll()fTT] are given by 



7^'(JV) = C A (4(N_ 2 - 2N_ - 2N 
7gg W = 4 Ca ( - 4 (N_ 2 - 2N_ - 2N+ + N +2 + 3) 



N +2 + 3) Si - y 



(3.1) 



Si -2 + Si, 2 + S 2 , 



+ ^(N+ - N +2 )S 2 - 4(N_ - 3N+ + N +2 + 1) 
+ y (N_ - N+)S 2 - I + 2Ss - y Si + 2S 3 



3S 2 - S 3 



109 /1VT 1VT 
+ 1 g-(N-+N+)Si 



(3.2) 



The three-loop gluon-gluon anomalous dimension reads, for n/ = 

■73091 



y$(N) = 16Ci((N_ 2 -2N_-2N + +N +2 +3) 



_5i-16Sl,-4 + ySl,-3 
+ 16Sl,-3,l + ^Sl.-2 + 451,-2,-2 - llSl,-2,l + 4Sl.-2,2 - T^Sl.l + 24Sl,l.-3 

6 108 

67 55 

+ HSl,l,-2 — I6S1, 1,-2,1 + 8Sl,l,3 g"^ 1 ' 2 8Sl,2,-2 + 8Sl,2,2 + y Si, 3 + 8Sl,3,l 



8S1, 
22 



395 
27 



11 



67, 



S 2 - 14S 2 ,-3 - — S 2 ,-2 + 8S 2 ,-2.1 - — S2.1 + 452,1,-3 + 852,1 



+ y 52,2 + 8S2,2,1 — 10S2,3 + 853,1,1 — 8S3,2 



9 

+ (2-N_ -N+) 



713 26 

— Si + ySl,_ 3 



145l,-2,l + Tr-51,-2 + — Sl.l - 145l,l,-2 + -T5-Sl,2 - 4Sl, 3 



+ (N_ - N-| 



473 
216 62 



12S 2 ,-3 + 552,-2 — 252,1 — 85 2 1 -2 + — S2.2 — 10S2,3 + -^7tS3 — 2053.-2 + — S3 1 

3 36 3 



16S 3 , 2 - 2IS4 - 26S 4 ,i 



(N_ - 



— — S2 + 8S2.-3 5"52,-2 — 8S 2 ,_2,1 

lOo 3 



8S 2 ,i,- 2 - ^5 2 ,2 - i-^S 3 + 8S 3 ,- 2 - ^S 3 ,i + 4S 3 , 2 + 44S 4 + 165 4 ,i - 8S 5 



+ (1-N 
1579 



52 — 8S2,-3 H — S 2 ,-2 + 852,-3,1 H — S2,l + 852,1,-2 + — 52,2 

10<5 3 9 3 



71 



— S3 + 8S 3 ,- 2 + — 5 3 ,i - 4S 3 , 2 - —Si - 16S 4 ,i + 36S 5 
lo o 6 



79 

32 



+ 4S_5 - 8S_4,i 
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+ ~~gS-3 — 45'_3 I -2 — 25-3,2 — 45*_2,-3 ^-S'_2,-2 + 45-2, -2,1 + 45'_2,l,-2 

16619 88 523 413 67 

" -igT Si ~ T Sl -~ 3 ~ Ts" 51 -- 2 + 115l -- a - 1 + To8 Sia _ 11Sl ' 1 -- 2 ~ Y Sl ' 2 

33 781 11 „ 67 „ 22 

— <~>1,3 + -Z-r"2 — 402.-3 + — 02,-2 + 4,52,-2,1 JT^.l + 4i>2,l,-2 77^2,2 

2 54 3 9 3 

+ ?5* 3 - 453,-2 + ¥-S a ,i - 2S 3 .2 - 85 4 ,i + 45 £ 
9 6 



(3.3) 



Note that harmonic sums up to weight 2n + 1 occur at order a™ +1 (N"LO). 
We stress that Eqs. and (ESJ 

are directly applicable only for even positive 
values of iV, while the lowest-order expression (|3.1[) holds for any positive integer. 
For general (non-integer) N, j^^N) can be obtained by numerically evaluating 

-l 



7 («)(7V) 



dx x N - 1 P {n \x) 



(3.4) 



using the a?-space results to which we turn now. 



3.2 Expressions in Bjorken-x space 

There is a theorem ensuring that the splitting functions P^ (x) can be uniquely 
reconstructed from their even-iV (or odd-./V) moments obtained in our calculations. 
In fact, the close relation between the harmonic sums and the harmonic polyloga- 
rithms facilitates an algebraic procedure |38I39| for the inverse Mellin transform. 

For a compact representation of the gluon-gluon splitting functions we use 
p gg {x) = (1 -x)' 1 +X- 1 -2 + x-x 2 
and an abbreviation for the harmonic polylogarithms |38j . 



H 



±(m+l),±(n+l) 



■Ho,...,o,±i,o,...,o,±i,...( x ) 



The one- and two-loop results |5l9j for = can then be written as 

11 



Pg\x) = C A (4 Pgg (x) + -5(1 
= 4C7 A 2 (2p gg ( a; ) 



-Ho,o 



67 
18 

44 



C2 + H ,o + 2Hi, + 2H ; 



2p g 



(3.5) 
C 2 + 2H_i,o 



U H + 8H 27 
— Mo + 8±lo,o — 



+ 27 - 12H 



+8(1 -x) 



+ 3Cs 



(3.6) 



The corresponding three-loop contribution reads 
P B W W = 16C A 3 ( Pgg (z) 



245 67 . 11 . 3.9 ., T „ TT . .__ 

— - —C.2 + y C3 - y^2 - 4H-3.0 + 6H_ 2 C2 + 4H_2,-1,0 

+— H_2,o — 4H_2,o.o — 4H_2,2 + —Ho — 7HoC3 + — Ho.o — 8Ho,o£2 + 4Ho,o,o,o — 6Hi^3 
3 6 9 

— 4Hi _2,o H — 7t- Hi ,o — 6Hi,o^ 2 + — Hi o.o + 8Hi ,o,o,o + 8Hno,o + 8Hi 2 ,o + 8Hi ,3 
9 6 
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+ -^H 2 - 4H2C2 + 10H 2 , ,o + 8H 2 ,i,o + 8H 2 , 2 + ^H 3 + 10H 3 ,o + 8H 3 ,i + 8H4 

9 D 



+Pgs(- X ) 



67 11 2 

y C2 + yC2 



411-3,0 + 16H-2C2 + 8H_2,-i,o — 18H_2,o,o — 12H_2, 2 

134, 



+I2H-1C3 + 8H_i,_2,o - 16H-i,_iC2 + 24H-i,_i,o,o + I6H-1 -1,2 - yE-y 

+18H-i,oC 2 - 16H_i,o,o,o - 4H_i,2,o - 16H_i, 3 - ^-HoC 2 - 5HoC 3 + ^H ,o - 8H ,oC 2 

6 9 



+4Hq, 0,0,0 + 2H2C2 + 2H3,o + 8H4 



+ 



(I-") 



16619 55 



11 TT , 67 TT 33 TT 67 TT 22 TT 

— H1C2 - yHi, + yHi, ,o - — H 2 + yH 2 ,o 



162 
+ 11 



C3 



H H t 413 H 
y HoC2 - ^Hi 



1 JL * 
— +X 

X 



389 2 
198 C2 " 3 H " 2 ' 



--H_iC^2 + H_i -1,0 — tt^tH-i + ttH-i,o,o + H_i,2 + —Ho — t:H 3 
2 198 3 54 6 



85 

Y C2 



,„„ , 6409 H +^H^ 1249 H 44H 44 H 110 TT ' 
+ooli_2,0 + ino Jlo + — — rlo,o — 44110,0,0 —^2,0 — H3 

lUo lo 00 



+(1 - X 
29 



11317 
108 
31, 



4H_ 3 ,o — 4H_2(^2 — 8H_2,-i,o 
25, 



19 

yH-2,0 



3 27 

3 Ho.o.o + ^Hi - -HxCa + yH 1;0 - yHi,o,o 



+ 1- 



- 12H_2,o,o y-Ho,o 

329, 11,, 

— C a + T (i + *)C8 



— — C2 2 — t-H_i<^2 — 3H_i -1.0 — H-i + 38H_i 0,0 + 25H_i 2 + — - Ho — 8H0C3 

5 2 6 216 

27 158 29 43 

+— H0C2 - 22H ,oC2 - — H 2 - 4H2C2 + — H 2 , + 10H 2 ,o,o - — H 3 + I6H3.0 + 26H 4 
2 9 3 d 

+ y C2 + 24Cs + 2C2 2 - 16H_3,o + 27H_ 2 ,o + y^H + y H0C2 - I6H0C3 - 29H ,o 



-4Ho,oC2 



40 



H ,o,o + 28rfi ,o,o,o + 27H 2 - 24H 2 ,o - 20H 3 + 5{1 - x) 



79 1, 

32 + 6 C2 



+ yC3 + ^<2 2 - 5Cs - C2C3 



(3.7) 



where, as in Eqs. I|3.5|) and (|3.6|) . all divergences for x — > 1 are to be read as +-dis- 
tributions. Functions H^ l (x) up to weight (number of indices) 2n occur at N"LO. 

A Fortran program for the harmonic polylogarithms up to weight four is 
available 001 • Nevertheless it is useful to have also more compact, if approximate 
representations of the three-loop splitting functions. Making use of the end-point 
behaviour discussed in the next section, Eq. I|3.7[) can be parametrized as [2] 

(x) = + 2643.521 V + 4425.894 5{1 - x) + 3589 L x - 20852 + 3968 x 
- 3363 x 2 + 4848 x 3 + L Li (7305 + 8757 L ) + 274.4 L Q 



where 



7471 Lo + 72 - 144 L% + 14214 x _1 + 2675.8 x ~ l L 



T>q — 1/(1 — x)+ , L\ — ln(l — x) , Lq = lnx 



(3-8) 



This parametrization deviates from the exact expression (|3.7(l by less than 0.1%, 
which should be perfectly sufficient for numerical applications. Note that the Mellin 
transform of Eq. (|3.8|) can be readily continued to complex values of N as required 
for the moment-space approach to the analysis of hard processes |41 | -|43 | . 
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4 Results and surprises for x — > 1 and x — > 

The end-point behaviour of the splitting functions is of particular interest. The 
leading contributions for x — > 1 are related to the cusp anomalous dimensions and 
are thus relevant beyond the context of parton distributions. The perturbativc 
stability at very small x, where potentially large ln fe x corrections occur, represents 
a much discussed topic directly relevant to analyses of collider processes. 



4-1 The large-x behaviour 

Up to N ra=2 LO, at least, the diagonal MS-scheme splitting functions are given by 



A 



n+1 



(1 



+ B„ a +1 5(1 - x) + C: +1 ln(l -x) + 0(1) . (4.1) 



In fact, the simple +-distribution 1/(1 — x) + constitutes the leading term to all 
orders — in contrast to, for example, the coefficient functions in DIS which include 
terms [(1 — a;) _1 ln fc (l — x)} + — and its coefficients A^ n form the perturbative 
expansion of the cusp anomalous dimensions of the respective Wilson lines |44| . 
For the quark case the known coefficients read, in our normalization, 



A? 



AC F 
8C F 



67 
18 



C 2 \Ca 



A? = 16 C F C; 



245 67 



11 



11 



!4~~¥ C2 + T C3 + T C2 





55 „ i. 


+ 16 C F n f 


V 2 C 3 

24 S3 



16 C F C A n. 



f 



209 
108 



10 , 



16 C F nf 



1 

27 



(4.2) 



The coefficients A^ are obtained from Eq. 1)4.2(1 by multiplication with Ca/C f . 
Note that the n^-independent parts of A^ consist of only one (the maximally non- 
abelian) colour factor as also predicted in ref. 03] . The = part of A^ and the 
complete Af are new results of refs. P and 0, respectively. Two computations of 
the n^-contribution to A£ were performed two years ago |20I22| . and its n^-part 
was already obtained in ref. ^3] . A^ in Eq. I|4.2|) agrees with the previous numerical 
estimate |45| which has been widely used in soft-gluon resummation analyses. 

It is also interesting that at three loops, as in the previous order, only a single 
logarithm occurs in Eq. 1)4.1)1 . Furthermore it turns out that there is an unexpected 
relation between the corresponding coefficients C^ +l and the A*, viz |1I2| 

C? = 0, C 2 a = (Af ) a , C* 3 a = 2A?A% . (4.3) 

The logarithmic structure of Eq. (|4.1() and especially the relation 1)4.3)1 call for an 
explanation and, possibly, a higher-order generalization. 

The large- a; limit of the quark-gluon and gluon-quark splitting functions reads 



,(«) 



V 2 ™" 1 D ab In 2 "" 



0(1) 



(4.4) 



See ref. [2] for the three-loop coefficients D, ah 
vanish for the choice Ca = C F 



2 i . Except for these coefficients 
rij: leading to a Af = 1 supersymmetric theory. 
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4-2 The small-x behaviour 

We start with the flavour non-singlet contributions which, in the present unpolar- 
ized case, are practically far less important then the singlet parts. However, the 
non-singlet small-x expansion includes two additional powers of In a; per order, i.e., 
terms up to ln 2 "x occur at N"LO. Thus the three- loop splitting functions 

P^ ±,s (x) = D^ln 4 x + ... + Dilnx + 0(1) (4.5) 

form a presently unique theoretical laboratory for studying the relative size of as 
many as four small-x logarithms as obtained from a complete (all-x) calculation. 

The numerical QCD values of the coefficients for i = — in Eq. 1|4.5[) read 

Dq Si 1.4321 

= 35.556- 3.1605 x^ 

D~ = 399.21 - 39.704 n f + 0.5926 nf 

Dg = 1465.9- 172.69^+4.3457x1/ . (4.6) 

See ref. for the analytic expressions and the similar case i — + . In both cases 
the leading coefficients Dq have been correctly predicted JB] on the basis of the 
resummation in ref. Note that the expansion 1)4.5(1 alternates for i = ±, and 
that the coefficients Dk rise sharply with k. For nj = 4 and i = — shown in Fig. 2, 

the modulus of the In 4 x (In 3 x) contribution is twice as large as that of the next 
term, In 3 x (In 2 x), only at extremely small x-values, x < 10 -14 (3 ■ 10 -10 ). 

The numerical situation is rather different for the d a bcd a bc contribution i = sP, 

£> S Si +1.4815 x^ , D{ Si -2.9630 x^ 

D| = +6.8918 x^ , D 3 S = +178.03 x^ . (4.7) 

Here the leading small-x terms do indeed provide a reasonable approximation, see 
Fig. 2. Note that the existence of a leading (ln 4 x) d a bcd a bc contribution for i = s 
is rather surprising. In fact, the presence of a leading small-x logarithm in a term 
unpredictable from lower-order structures appears to call into question the very 
concept of the small-x resummation of the double logarithms a™ +1 In 2 " x. 

The leading small-x terms of the three-loop singlet splitting functions are 

P&l^oto = E i h ^ + E 2 h ~ + °(^ X ) ■ ( 4 - 8 ) 

In general, x -1 ln fc x contributions with k < n occur in P gq and P gg at N ra LO. The 
highest of these terms (k = n) have, however, vanishing coefficients for n= 1, 2 as 
predicted by the leading-logarithmic BFKL equation [47148] . see also ref. |49| . 
In QCD the numerical values of the coefficients in Eq. I|4.8|l are given by 

E m ^ -132.74^ , £ 2 qq Si -506.00x^ + 3.1605x1/ 

E? B Si -298.67^ , E. qs = -1268.3x^ + 4.5761x2,/ 

E? q Si 1189.3+ 71.083 n f , Ef q Si 6163.1-46.408x^-2.3704x2./ 

Ef s = 2675.9+ 157.27 n f , Ef g Si 14214. + 182.96 n f - 2.7984 n/ . (4.9) 
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The analytical results can be found in ref. 0. The coefficients i? 1 qa and E^ g agree 
with those derived from the small- a; rcsummation in refs. |17j and |19j . respectively, 
after transforming the latter result to the MS scheme jHOJ- Ef^ was unknown before 
ref. 0. For n f = 3 ... 5 the ratios E^/Ef 3 are 3.7. . .4.7. Thus the corrections 
due to the non-logarithmic 1/x terms amount to less than 50% only at x < 10 -4 . 



1 1 


1 1 1 


\ 


(l-x)Pl 2) (x) 


\ 




- \ 
' \ 




\ -\ 

\\ • V 
\\ ■ 

N >v ■ . 




s \. 


1 1 1 





p s (2) (x) ; 






-\ (X 0.25) 


exact 

• ■ ■ ln 4 x 




+ln 3 x : 

+ In 2 x - 

- - + In 1 x I 








...i i " 



3000 




Figure 2. The small- £ behaviour of non-singlet (upper row) and singlet NNLO splitting functions 
for rij. =4. Also shown are the (successive) approximations (14.51 and 14.81 by the leading small- a; 
logarithms and, for the singlet cases, the errors bands l"i 1 1 used in the analyses of refs. [52 53 . 
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5 The size of the corrections 

Finally we discuss the numerical effects of our new contributions to Eq. 11.31) . 
For brevity, we confine ourselves to four quark flavours and a typical scale u 2 ~ 
30 ... 50 GeV 2 . We employ an order- independent value of the strong coupling, 

a s (/4 n f = 4) = 0.2 , (5.1) 

facilitating a direct comparison of the MS evolution kernels at LO, NLO and NNLO. 

5.1 N-space: anomalous dimensions 

The singlet anomalous dimensions 7 ff , are displayed in Fig. 3 for the standard 
choice \i r = fif of the renormalization scale tacitly made already in sections 3 
and 4. The NNLO corrections are markedly smaller than the NLO contributions. 
For the choice l|5.1|) they amount, at N>2, to less than 2% and 1% for the large 
diagonal quantities 7 qq and 7 gg , respectively, while for the much smaller off-diagonal 
anomalous dimensions 7 qg and 7 gq values of up to 6% and 4% are reached. Also 
shown in Fig. 3 is the pure-singlet contribution 7 ps defined by Eq. I|2.8|l . At N > 2 
this quantity receives very large relative (but tiny absolute) NNLO corrections. 




N N 

Figure 3. The perturbativc expansion of the singlet anomalous dimensions 7 ff , up to NNLO. 
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5.2 Scale derivatives of x- space parton distributions 

In Figs. 4 and 5 wc show the logarithmic derivatives fk = din fk/dln fxj for the 
sufficiently realistic - and like Eq. (|5.1|l order-independent - model distributions 

xq s (x,fx£) = 0.6x- a3 (l-x) 35 (l + 5.0a; a8 ) 

xg(x,^) = 1.6x- a3 (l-x) 45 (l-0.6a; ' 3 ) . (5.2) 

At large x the NNLO corrections to the non-singlet evolution illustrated in Fig. 4 
are very similar for all three cases (12. 5[) . They amount to 2% or less for x > 0.2. 
thus being smaller than the NLO contributions by a factor of about eight. For 
q~ s the same suppression is also found in the region 10~ 5 < x < 10~ 2 . The NNLO 
effects are even smaller for g + (not shown) at small x, but considerably larger for 
q^s at x < 10~ 3 due to the additional effect of the new quantity P n s s in Eq. (|2.6|l . 




X X 

Figure 4. The perturbative expansion up to NNLO of the factorization-scale derivatives q na ' v for 
the initial conditions 15.11 . 15.21 and the standard choice fj, r = /if of the renormalization scale. 

Also for the singlet quark distribution (upper row of Fig. 5) the ratio of the NLO 
and NNLO corrections is about eight over the full x-rangc. However, at small x 
- where P qg * g dominates in Eq. I|2.7|l - the LO results are anomalously small as 
P qq and P qg , unlike at higher orders, do not include 1/x terms at first order. The 
situation is quite different for the evolution of the gluon density (dominated by 
Pgg * g at all x). Here the NLO contributions appear atypically small at low x, 
cf. the remark below Eq. I|4.8() . Thus the ratio of the NNLO and NLO corrections 
is rather large here, despite the former amounting to only 3% for x as low as 10 ~ 4 . 

It is also interesting to consider the stability of the above results under variations 
of the renormalization scale fi r . For pi r ^ fj,f the perturbative expansion of the 
splitting functions up to NNLO reads, with Lf r = ln(/ij//x 2 ), 

P = a s (^ 2 )p(°) + a*(rf)(pU-f3 pML fr ) 

+ a?(rf) (P^ {/?1 p(°) + 2/? P^}L fr + 01 P<®L) r ) . (5.3) 
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d In q s / d In J 




d In g / d In \xj 




0.9 
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1 1 1 
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7 NLO/LO 

! NNLO/NLO 


■ 
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X 



10 



Figure 5. As Fig. 4, but for the singlet quark and gluon distributions q B and g. The spikes close 
x — 0.1 in the right plots reflect the zeros of /j, and do not constitute large absolute corrections. 

Here a s (//j5) is obtained at N™LO from the value l|5.1|) at the scale by solving 



din ^ 



/?(a») 



1=0 



(5.4) 



The MS expansion coefficients (3k are presently known up to k = 3 |54j-|58j. 

The dependence of the above results on fi r can be found in Fig. 8 of ref. pQ and 
Figs. 9 and 10 of ref. [5] for selected values of x. Here we only show, in Fig. 6, 
the relative scale uncertainties A/fc of the average /^-derivatives as conventionally 
estimated by varying \i r up to a factor of two with respect to/if, 



A/ = 



max / — min / 
2 | average /| 



(5.5) 



For the valence, (singlet-quark, gluon) distributions, these uncertainty estimates 
amount to 3% (3%, 1%) or less at x > 5 • 10~ 4 (5 • 10~ 3 , 3 • 10~ 4 ), an improvement 
by more than a factor of three with respect to the corresponding NLO results. 
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0.15 



0.05 




0.15 



0.(15 



Figure 6. The renormalization scale uncertainties of the NLO and NNLO predictions for the evo- 
lution of the non-singlet and singlet distributions as estimated by the quantities A/fe in Eq. 15.51 . 



6 Summary and outlook 

We have calculated the complete third-order splitting functions for the evolution 
of unpolarized parton distributions in perturbativc QCD. The computation is per- 
formed in Mellin-iV space and follows the previous fixed- N computations |12|-|14| 
inasmuch as we compute the partonic structure functions in deep-inelastic scattering 
at even or odd N. Our calculation, however, provides the complete TV-dependence 
from which the splitting functions in Bjorken-x space can be uniquely reconstructed. 

A salient feature of our approach is that it facilitates very efficient checks of 
our elaborate schemes for the reduction of the required three-loop integrals by 
using the Mincer program |24I25| . By keeping terms of order e° in dimensional 
rcgularization throughout the calculation, we have also obtained the third-order 
coefficient functions for the structure functions F2 and Fl in electromagnetic and 
for F 3 in charged-current DIS [231 - The present method can be used to generalize 
our fixed- TV three-loop calculation of the photon structure [HH] to all N. It should 
also be possible to obtain the polarized NNLO splitting function in this manner. 
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Our results completely agree with all partial results available in the literature 
for fixed moments large- limits |15I16| . small- a; behaviour E7j-^5j and 

large- a; structure |22I44| . With the (relatively unimportant) exception of P gq shown 
above and P n s s they also fully agree with the uncertainty bands of ref. |SJ used in 
provisional NNLO analyses 52 53 . Those analyses thus remain valid. The results 
do, however, exhibit some unexpected features, most notably a suggestive relation 
between large- a; coefficients and the presence of a leading small- a; term in the new 
dabedabe contribution P n s s to the splitting function for the valence distribution. 

The effect of the three-loop (NNLO) corrections on the evolution of the parton 
densities is small at x > 10~ 3 . For a s = 0.2, for example, both the corrections and 
the \x r variation amount to less than 2% at large x; and the NNLO effects are 
about eight times smaller than the NLO contributions, implying that the evolution 
is perturbatively stable down to rather low scales. For x < 10~ 3 the corrections 
increase with decreasing x. As the knowledge of the leading small- a; terms is not 
sufficient, further improvements in this region would require considerable efforts, 
including at least an extension of the Mincer program [2 to four loops. 
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